Ce document montre l'essentiel des manipulations Ơ effectuer avec R pour rƩpondre aux questions du TP1 d'Analyse UnivariƩe. Il suppose une connaissance prƩalable de la syntaxe R, de prƩfƩrence Ơ l'aide de l'IDE RStudio et les bases du package sf (tutoriel ici) qui implƩmente la norme SF comme dans PostGIS.

La connaissance du package dplyr et de sa notation Ơ base de l'opƩrateur pipe (%>%) sera utile pour Ʃcrire des chaƮnes de traitements plus concises que celles listƩes ici.

0.1 Chargement des donnƩes.

Nous avons rƩcupƩrƩ le fichier des donnƩes spatialisƩ sous le format GeoJSON. disponible ici.

D'autres formats sont disponibles, notamment via l'API du site. Nous allons utiliser le package sf pour manipuler les donnƩes spatiales. Pour lire le GeoJSON , nous avons besoin d'installer un package spƩcifique : geojsonsf, qui se charge de la conversion entre des objets geographiques au format GeoJSON et leur reprƩsentation au format sf.

0.1.1 Installation des librairies

L'installation et le chargement des packages au moyen des fonctions install.packages("nom_du_package") et library(nom_du_package)

install.packages("geojsonsf")
library(geojsonsf)
install.packages("sf")
library(sf)

Nous installons et chargeons aussi les packages dplyr (manipulation aisƩe de donnƩes) et cartography (cartographie correcte des objets spatiaux).

install.packages("dplyr")
library(dplyr)
install.packages("cartography")
library(cartography)

0.1.2 Lecture de fichier

Le fichier geoJSON est chargƩ Ơ l'aide de la fonction geojson_sf().


N.B. la fonction geojson_sf requiert le chemin d'accès au fichier, fourni sous sa forme relative ou absolue. Le repertoire de travail de R est défini par la commande setwd("/chemin/vers/le/fichier"). Pour connaître le repertoire de travail courant, utiliser la fonction getwd()

N.B. le nombre de lignes et les valeurs de données qui apparaissent dans ce support peuvent varier, les données étant mises à jour régulièrement sur le site opendata.paris.fr


arbres <- geojson_sf("les-arbres.geojson")
names(arbres) # affiche le nom des colonnes (variables)
##  [1] "idbase"             "remarquable"        "genre"             
##  [4] "stadedeveloppement" "arrondissement"     "idemplacement"     
##  [7] "geo_point_2d"       "geometry"           "adresse"           
## [10] "libellefrancais"    "complementadresse"  "domanialite"       
## [13] "circonferenceencm"  "typeemplacement"    "hauteurenm"        
## [16] "varieteoucultivar"  "espece"

Le chargement du fichier peut être un peu long: il contient 203818 observations le jour de la rédaction de ce support. Pour connaitre la structure de l'objet , nous utilisons la fonction str() qui nous donne un aperçu du type des colonnes du dataframe que nous manipulons. C'est un dataframe particulier puisqu'il est également de la classe sf : l'un de ses colonnes est dédié à la description de sa géométrie (ici, des points 2D) dans la colonne dédiée : geometry

0.2 AperƧu des donnƩes

La fonction head() permet d'observer les \(n\) premières lignes d'un dataframe. Ici la fonction paged_table est utilisée pour un affichage dynamique plus pratique du tableau dans ce support HTML.

paged_table(head(arbres, 10))

0.3 Projeter et Afficher les donnƩes

La librairie sf surcharge la fonction d'affichage par défaut de R , plot, pour afficher la géométrie des données comme un objet géographique et non comme un nuage de points ou une courbe. Auparavant , nous devons fixer le système de coordonnées de références de l'objet, pour que les données soient correctement projetées à l'affichage.

La fonction st_crs() appliquƩe sur un objet spatial retourne son CRS s'il est dƩfini , ou permet de fixer avec l'opƩrateur d'affectation <-. Nous allons utiliser Lambert 93 (EPSG 2154). L'opƩration se fait en deux temps : transformation selon le CRS puis affectation de l'attribut.

st_crs(arbres) 
## Coordinate Reference System:
##   User input: 4326 
##   wkt:
## GEOGCS["WGS 84",
##       DATUM["WGS_1984",
##         SPHEROID["WGS 84",6378137,298.257223563,
##           AUTHORITY["EPSG","7030"]],
##         AUTHORITY["EPSG","6326"]],
##       PRIMEM["Greenwich",0,
##         AUTHORITY["EPSG","8901"]],
##       UNIT["degree",0.0174532925199433,
##         AUTHORITY["EPSG","9122"]],
##       AXIS["Latitude",NORTH],
##       AXIS["Longitude",EAST],
##     AUTHORITY["EPSG","4326"]]
arbres <- st_transform(arbres,2154)
arbres <- st_set_crs(arbres, 2154)
st_crs(arbres)  # nouveau CRS de l'objet aprĆØs transormation
## Coordinate Reference System:
##   User input: EPSG:2154 
##   wkt:
## PROJCS["RGF93 / Lambert-93",
##     GEOGCS["RGF93",
##         DATUM["Reseau_Geodesique_Francais_1993",
##             SPHEROID["GRS 1980",6378137,298.257222101,
##                 AUTHORITY["EPSG","7019"]],
##             TOWGS84[0,0,0,0,0,0,0],
##             AUTHORITY["EPSG","6171"]],
##         PRIMEM["Greenwich",0,
##             AUTHORITY["EPSG","8901"]],
##         UNIT["degree",0.0174532925199433,
##             AUTHORITY["EPSG","9122"]],
##         AUTHORITY["EPSG","4171"]],
##     PROJECTION["Lambert_Conformal_Conic_2SP"],
##     PARAMETER["standard_parallel_1",49],
##     PARAMETER["standard_parallel_2",44],
##     PARAMETER["latitude_of_origin",46.5],
##     PARAMETER["central_meridian",3],
##     PARAMETER["false_easting",700000],
##     PARAMETER["false_northing",6600000],
##     UNIT["metre",1,
##         AUTHORITY["EPSG","9001"]],
##     AXIS["X",EAST],
##     AXIS["Y",NORTH],
##     AUTHORITY["EPSG","2154"]]

Nous pouvons maintenant afficher les donnƩes avec la fonction plot

plot(arbres)
## Warning: plotting the first 9 out of 16 attributes; use max.plot = 16 to plot
## all

C'est très long ! Celà est dû au fait que la fonction plot affiche autant de cartes que de variables (colonnes) de l'objet spatial, jusqu'à un maximum de 9 colonnes apr défaut.

Pour le moment nous n'allons que tracer la gƩomƩtrie des donnƩes : l'attribut geometry de l'objet

plot(arbres$geometry, cex = 0.01)

Tout autre attribut peut être choisi , la fonction de dessin de R se charge de trouver une échelle de couleur en fonction des valeurs que prend la varaible.

plot(arbres["arrondissement"], cex = 0.01, graticule=T)

La sequence de commande suivantes transforme l'attributs arrondissement en facteurs et filtre les donnƩes de faƧon Ơ ne conserver que les arrondissements de Paris intra muros.

Nous utilisons les fonctions:

arbres$arrondissement <- as.factor(arbres$arrondissement)
levels(arbres$arrondissement) #  modalitƩs de la variable
##  [1] "BOIS DE BOULOGNE"  "BOIS DE VINCENNES" "HAUTS-DE-SEINE"   
##  [4] "PARIS 10E ARRDT"   "PARIS 11E ARRDT"   "PARIS 12E ARRDT"  
##  [7] "PARIS 13E ARRDT"   "PARIS 14E ARRDT"   "PARIS 15E ARRDT"  
## [10] "PARIS 16E ARRDT"   "PARIS 17E ARRDT"   "PARIS 18E ARRDT"  
## [13] "PARIS 19E ARRDT"   "PARIS 1ER ARRDT"   "PARIS 20E ARRDT"  
## [16] "PARIS 2E ARRDT"    "PARIS 3E ARRDT"    "PARIS 4E ARRDT"   
## [19] "PARIS 5E ARRDT"    "PARIS 6E ARRDT"    "PARIS 7E ARRDT"   
## [22] "PARIS 8E ARRDT"    "PARIS 9E ARRDT"    "SEINE-SAINT-DENIS"
## [25] "VAL-DE-MARNE"
arbres_intramuros <-  filter(arbres, !(arrondissement %in% c("BOIS DE BOULOGNE",  "BOIS DE VINCENNES", "HAUTS-DE-SEINE", "SEINE-SAINT-DENIS", "VAL-DE-MARNE")))
# on retire les valeurs modales inusitƩes de la variable 
arbres_intramuros$arrondissement <- as.character(arbres_intramuros$arrondissement)
arbres_intramuros$arrondissement <- as.factor(arbres_intramuros$arrondissement)
levels(arbres_intramuros$arrondissement)
##  [1] "PARIS 10E ARRDT" "PARIS 11E ARRDT" "PARIS 12E ARRDT" "PARIS 13E ARRDT"
##  [5] "PARIS 14E ARRDT" "PARIS 15E ARRDT" "PARIS 16E ARRDT" "PARIS 17E ARRDT"
##  [9] "PARIS 18E ARRDT" "PARIS 19E ARRDT" "PARIS 1ER ARRDT" "PARIS 20E ARRDT"
## [13] "PARIS 2E ARRDT"  "PARIS 3E ARRDT"  "PARIS 4E ARRDT"  "PARIS 5E ARRDT" 
## [17] "PARIS 6E ARRDT"  "PARIS 7E ARRDT"  "PARIS 8E ARRDT"  "PARIS 9E ARRDT"

On peut maintenant tracer les positions des arbres de Paris intra-muros, après avoir vérifié le CRS de cette couche par la fonction st_crs() (pour superposer plus tard les deux couches si celles-ci )

st_crs(arbres_intramuros) 
## Coordinate Reference System:
##   User input: EPSG:2154 
##   wkt:
## PROJCS["RGF93 / Lambert-93",
##     GEOGCS["RGF93",
##         DATUM["Reseau_Geodesique_Francais_1993",
##             SPHEROID["GRS 1980",6378137,298.257222101,
##                 AUTHORITY["EPSG","7019"]],
##             TOWGS84[0,0,0,0,0,0,0],
##             AUTHORITY["EPSG","6171"]],
##         PRIMEM["Greenwich",0,
##             AUTHORITY["EPSG","8901"]],
##         UNIT["degree",0.0174532925199433,
##             AUTHORITY["EPSG","9122"]],
##         AUTHORITY["EPSG","4171"]],
##     PROJECTION["Lambert_Conformal_Conic_2SP"],
##     PARAMETER["standard_parallel_1",49],
##     PARAMETER["standard_parallel_2",44],
##     PARAMETER["latitude_of_origin",46.5],
##     PARAMETER["central_meridian",3],
##     PARAMETER["false_easting",700000],
##     PARAMETER["false_northing",6600000],
##     UNIT["metre",1,
##         AUTHORITY["EPSG","9001"]],
##     AXIS["X",EAST],
##     AXIS["Y",NORTH],
##     AUTHORITY["EPSG","2154"]]
plot(arbres_intramuros["arrondissement"], cex=0.1)

0.4 DonnƩes vectorielles des quartiers de Paris.

Nous allons maintenant charger les contours des quartiers de Paris, disponibles sur (https://opendata.paris.fr/explore/dataset/quartier_paris/information/). Chaque arrondissement contient 4 quartiers, on a donc 80 unitƩs spatiales Ơ considƩrer.

quartiers <-(read_sf("quartier_paris.shp"))
quartiers <- st_transform(quartiers, 2154)
quartiers <-  st_set_crs(quartiers,2154)

plot(quartiers$geometry)
plot(arbres_intramuros["arrondissement"], add=TRUE, alpha=0.5, cex=0.1  )

0.5 "Cartographie" simple

Pour tracer un objet spatial , il suffit d'appliquer la fonction plot() sur cet objet. Cette fonction peut Ʃgalement colorer la gƩomƩtrie de l'objet en fonction d'une variable.

0.5.1 Les 6 genres d'arbres les plus reprƩsentƩs

Nous utilisons les fonctions count()du package dplyr, pour compter les nombre d'arbre par genre. (c'est l'Ʃquivalent d'un group_by, suivi d'une aggregation de comptage : count) Nous utilisons la fonction top_n du package dplyr pour selectionner les 6 genres les plus reprƩsentƩs. Nous utilisons la fonction filter et l'operateur %in% pour ne conserver que les arbres de ces six genres.

count_by_genre <-  count(arbres_intramuros,genre, libellefrancais, sort = T,name = "nb_indiv")
head(count_by_genre, 10)
## Simple feature collection with 10 features and 3 fields
## geometry type:  MULTIPOINT
## dimension:      XY
## bbox:           xmin: 644413.8 ymin: 6857489 xmax: 657170.6 ymax: 6867050
## CRS:            EPSG:2154
## # A tibble: 10 x 4
##    genre   libellefrancais nb_indiv                                     geometry
##    <chr>   <chr>              <int>                             <MULTIPOINT [m]>
##  1 Platan… Platane            37758 ((645035 6861100), (645068.5 6861382), (645…
##  2 Aescul… Marronnier         20022 ((644413.8 6861116), (644420.5 6861114), (6…
##  3 Tilia   Tilleul            16918 ((645059.1 6860934), (645070.1 6861108), (6…
##  4 Acer    Erable             12455 ((645032.2 6860884), (645037.9 6860896), (6…
##  5 Sophora Sophora            10752 ((645069.1 6861372), (645088.8 6860152), (6…
##  6 Pinus   Pin                 3719 ((645037.9 6860889), (645057.3 6861095), (6…
##  7 Celtis  Micocoulier         3693 ((645102.9 6860998), (645167.5 6861102), (6…
##  8 Fraxin… FrĆŖne               3551 ((645068.1 6860577), (645074.4 6861408), (6…
##  9 Prunus  Cerisier Ć  fle…     3456 ((645037.9 6860884), (645056.4 6860576), (6…
## 10 Pyrus   Poirier Ć  fleu…     3281 ((645093.6 6861056), (645103.2 6861005), (6…
six_genres <- top_n(count_by_genre,6,nb_indiv)
head(six_genres)
## Simple feature collection with 6 features and 3 fields
## geometry type:  MULTIPOINT
## dimension:      XY
## bbox:           xmin: 644413.8 ymin: 6857489 xmax: 657126.1 ymax: 6867050
## CRS:            EPSG:2154
## # A tibble: 6 x 4
##   genre   libellefrancais nb_indiv                                      geometry
##   <chr>   <chr>              <int>                              <MULTIPOINT [m]>
## 1 Platan… Platane            37758 ((645035 6861100), (645068.5 6861382), (6450…
## 2 Aescul… Marronnier         20022 ((644413.8 6861116), (644420.5 6861114), (64…
## 3 Tilia   Tilleul            16918 ((645059.1 6860934), (645070.1 6861108), (64…
## 4 Acer    Erable             12455 ((645032.2 6860884), (645037.9 6860896), (64…
## 5 Sophora Sophora            10752 ((645069.1 6861372), (645088.8 6860152), (64…
## 6 Pinus   Pin                 3719 ((645037.9 6860889), (645057.3 6861095), (64…

L'immense avantage du package sf : les gƩomƩtries ont ƩtƩ conservƩes lors des manipulation de regrouppement , filtrage etc. Il n'est donc pas nƩcessaire de refaire les opƩrations sur les gƩomƩtries , qui ont "suivi" les rƩsultats lors des opƩrations count et top_n .

plot(six_genres["libellefrancais"], 
       key.pos=1, cex=0.2,
       key.length = 0.9, 
       main="Les six genres d'arbres les plus reprƩsentƩs Ơ Paris")

Que peut-on dire de l'implantation de ces 6 genres ?

0.5.2 La domanialitƩ

Afficher cette variable est immƩdiat. La lƩgende est perfectible pour le moment; pour faire une vraie carte , voir Ơ la fin du support l'utilisation du package cartography.

 plot(arbres_intramuros["domanialite"], 
      key.pos=4, cex=0.2, 
      key.length = 1,
      key.width = lcm(4.5),
      main="DomanialitƩ des arbres de Paris")

0.6 Calcul du nombre d'arbres par quartier

0.6.1 Aggregation spatiale

On rƩalise une jointure (spatiale) entre quartiers et arbres_intramuros, pour inclure les informations de quartiers aux arbres Ơ l'aide du prƩdicat spatial st_within

arbres_in_quartiers <- st_join(arbres_intramuros,quartiers, join=st_within)
nrow(arbres_in_quartiers) 
## [1] 164612

On obtient alors un nouveau dataframe spatial, contenant les 163367 arbres situƩs dans les quartiers intramuros, et augmentƩs des variables issues du dataframe quartier.

On peut alors compter le nombre d'arbres par quartiers avec la fonction table qui compte (entre autres) le nombre d'individus (ligne) par valeurs distinctes d'une de ses variables. On peut Ʃtendre le nombre de variables pour obtenir des tables de contingence.

nb_arbres_by_quartier <- table(arbres_in_quartiers$c_qu)

Pour ajouter cette information au dataframe quartiers, l'une des possibilitƩ est de trier le dataframe par son code de quartier c_qu puis d'ajouter la colonne nb_arbres_by_quartiers obtenue prƩcƩdemment.

on peut par exemple utiliser la fonction arrange du package dplyr

quartiers <-  arrange(quartiers,c_qu)
quartiers$nb_arbres <-  nb_arbres_by_quartier

Avant de cartographier cette valeur, regardons sa distribution à l'aide d'un histogramme pour déterminer si un traitemment particulier doit être envisagé (en cas de distribution particulièrement biscornue)

hist(quartiers$nb_arbres,breaks = 10)

0.7 Carte choroplĆØte du nombre d'arbres (mauvaise pratique)

plot(quartiers["nb_arbres"], main=NULL, key.pos = 4)
title("Nombre d'arbres par quartier administrif de Paris",sub = "ReprƩsentation sƩmiologiquement discutable" )

0.8 "Raster" du nombre d'arbres

Pour réaliser une sorte de raster, nous allons créer une grille qui couvre la zone d'étude, et projeter les points du semis (le dataframe arbres_intramuros) dans les cellules créées (c'est également une façon d'approcher une carte de densité 2D de façon discrète)

Afin de ne pas faire une grille qui couvre la boƮte englobante de la zone d'Ʃtude, mais que la grille ne couvre que les gƩomƩtries de chaque quartier, on rƩalise une intersection entre la grille raster et l'enveloppe des quartiers.

N.B. ce n'est pas la faƧon canonique de rƩaliser un raster, pour cela il faut utiliser la library raster et rƩaliser le raster Ơ partir l'aide de la fonction

rast1 <- st_make_grid(quartiers, square = TRUE, n=40, what="polygons") ## raster de 1600 cellules
plot(quartiers$geometry)
plot(rast1,  add=T)

#restriction du raster Ơ la geomƩtrie des quartiers par intersection avec son enveloppe
envelop_qu <-  st_union(quartiers)
plot(envelop_qu)

rast2 <-  st_intersection(envelop_qu, rast1) ## raster de 1065 cellules

# on retransforme l'objet rast2 qui est une collection (sfc) en objet sf
rast2 <- st_sf(rast2)
plot(rast2)

#on récupère les points contenus dans les cellules
predicat_intersect <- st_contains(rast2,arbres_intramuros)

#on compte le nombre de points par cellules avec la fonction length
rast2$nb_arbres <- unlist(lapply(predicat_intersect, length)) 
hist(rast2$nb_arbres)

plot(rast2)

A quoi correspondent les cellules jaunes ?

On peut rƩitƩrer le processus en changeant la taille de la grille ou le type de grille (hexagonale)

rast1 <- st_make_grid(quartiers, square = TRUE, n=80)
rast2 <-  st_intersection(envelop_qu, rast1) 
rast2 <- st_sf(rast2)
predicat_intersect <- st_contains(rast2,arbres_intramuros)
rast2$nb_arbres <- unlist(lapply(predicat_intersect, length)) 

plot(quartiers$geometry, border="white", bgc="#222222"  )
plot(rast2, border=NA, key.pos=4, add=T)
plot(quartiers$geometry, border="white", bgc="#222222", add=T, alpha=0.8,lwd=0.3,key.pos=4)

rasthex <- st_make_grid(quartiers, square = FALSE, n=60, what="polygons")
rasthex <- st_sf(rasthex)
#les hexagones s'arrètent à l'evloppe, pas besoin de la créer
predicat_intersect <- st_contains(rasthex,arbres_intramuros)
rasthex$nb_arbres <- unlist(lapply(predicat_intersect, length)) 

plot(quartiers$geometry,bgc="#222222")
plot(rasthex,  key.pos=4, add=T, lwd=0.2)
plot(quartiers$geometry, border="white", bgc="#222222", add=T, alpha=0.8,lwd=0.3,key.pos=4)

On peut noter que le changement de dƩcoupage de l'espace dƻ Ơ la taille ou la forme des cellules change l'aspect de la carte.

0.9 Calcul et cartographie de la densitƩ d'arbres par quartiers

0.9.1 Calcul et carte choroplète de la densité

Pour calculer la densitƩ d'arbres par quartier, nous devons compter le nombre d'arbres par quartier, et diviser ce nombre par la surface du qartier.

Nous avons dƩjƠ ajoutƩ un attribut nb_arbres Ơ l'objet quartiers. L'objet quartiers contient aussi un attribut surface, issu des donnƩes brutes. On doit d'abord vƩrifier que cette valeur est correcte (la projection des donnƩes est Ʃquivalente) en calculant l'aires des quartiers et en les comparant Ơ l'attribut prƩ-existant

ecarts <-  as.numeric(st_area(quartiers)) - quartiers$surface
ecarts
##  [1] 1.866755e-05 8.272764e-06 5.813083e-06 5.191891e-06 4.138594e-06
##  [6] 4.982692e-06 5.938637e-06 6.454240e-06 7.424154e-06 6.026763e-06
## [11] 7.710827e-06 3.849418e-06 6.659189e-06 9.832554e-06 1.031480e-05
## [16] 8.575618e-06 1.227891e-05 1.769396e-05 1.484947e-05 9.281968e-06
## [21] 6.717455e-06 1.572375e-05 1.857313e-05 6.059418e-06 1.743273e-05
## [26] 2.260786e-05 1.833355e-05 2.937787e-05 2.528704e-05 1.730002e-05
## [31] 1.635018e-05 2.579298e-05 1.576392e-05 1.172686e-05 9.100710e-06
## [36] 1.092535e-05 2.064311e-05 9.532028e-06 1.345109e-05 1.900073e-05
## [41] 1.540396e-05 1.834927e-05 2.513756e-05 2.110482e-05 1.282701e-04
## [46] 1.555709e-04 4.194141e-05 2.657715e-05 2.524280e-05 6.549200e-05
## [51] 4.865695e-05 1.541385e-05 2.415944e-05 2.837135e-05 2.953457e-05
## [56] 3.841263e-05 5.984306e-05 3.550434e-05 3.195251e-05 5.530147e-05
## [61] 1.364052e-04 1.210235e-04 6.643729e-05 3.027916e-05 3.236555e-05
## [66] 2.936390e-05 3.111642e-05 3.004400e-05 4.075211e-05 3.607967e-05
## [71] 2.370775e-05 2.878439e-05 2.853689e-05 5.197572e-05 3.966014e-05
## [76] 2.765586e-05 1.754111e-05 3.255368e-05 3.333972e-05 4.542875e-05

Les écarts sont très faibles (de l'ordre de \(10^{-5}m^2\)) , suffisamment pour que l'attribut surface soit directement utilisé pour le calcul de la densité.

On crƩe une variable dens dans l'objet quartier dont la valeur est le ratio entre la variable surface et la variable nb_arbres. On aura vƩrfiƩ auparavant qu'il n'y a pas de valeurs manquantes dans l'objet , avec la fonction anyNA()

anyNA(quartiers)
## [1] FALSE
quartiers$dens <- quartiers$nb_arbres / quartiers$surface
plot(quartiers["dens"], main="DensitƩ d'arbres par quartier")

Connaissant l'implantation des arbres et la formes des quartiers , que dire des densitƩ des quartiers sud-ouest et sud-est ?

Que nous dit la comparaison entre la carte choroplète de nombre et celle de densité ?

Quel est le quartier le plus dense en arbres ?

leplusdense <- top_n(quartiers,1,dens) 
leplusdense$l_qu
## [1] "PĆØre-Lachaise"

La densité des bois de Paris ne vous choque pas ? D'où viennent ces faibles densités ? Comment les corriger ?

0.9.2 Cartographie par symboles proportionnels

Pour utiliser des symboles proportionnels , nous devons utiliser un package spƩcifique : cartography. Dans ce package , l'ajout de couche au dessin courant est activƩ par dƩfaut (add=TRUE)

library(cartography)
plot(quartiers$geometry, bgc="#888888", border="white", lwd=0.3)
propSymbolsLayer(quartiers,var = "nb_arbres", inches = 0.15, legend.pos = NA)
# couches de disposition des lƩgendes et du texte
layoutLayer(title = "Nombre d'arbres Ć  Paris par quartier",  
            frame =F, north = TRUE, author="M2 IGAST2019-2020",
            sources = "Opendata.paris.fr 2019", 
            scale = 50)
legendCirclesSymbols(pos = "topleft", title.txt = "nombre d'arbres", title.cex = 0.8, cex = 1,
  border = "black", lwd = 1, values.cex = 0.6, var=c(50,1000,3500,7000), inches=0.15,
  col = "#E84923", frame = FALSE, values.rnd = 0, style = "e")

Le package cartography est un excellent package, que je recommande vivement.